///>£ -/-/I 


(AIAA 99-3177, AIAA 17th Applied Aerodynamics Conference, Norfolk, VA, June 1999) 

COMPUTATIONAL MODELING FOR THE FLOW OVER <0 
A MULTI ELEMENT AIRFOIL 


William W. Liou* and Fengjun Liu** 
Department of Mechanical and Aeronautical Engineering 
Western Michigan University 
Kalamazoo, MI 49008 


Abstract 

The flow over a multi-element airfoil is computed 
using two two-equation turbulence models. The com- 
putations are performed using the INS2D Navier- 
Stokes code for two angles of attack. Overset grids 
are used for the three-element airfoil. The computed 
results are compared with experimental data for the 
surface pressure, skin friction coefficient, and velocity 
magnitude. The computed surface quantities generally 
agree well with the measurement. The computed re- 
sults reveal the possible existence of a mixing-layer-like 
region of flow next to the suction surface of the slat for 
both angles of attack. 

I. Introduction 

Flows past multi-element airfoils have been sub- 
jected to intensive experimental and computational 
studies 1,2,3 for the past two decades. An accurate pre- 
diction of such flows can enhance the performance and 
the safety factor of aircrafts in high-lift operations. 
Even with the advances of computational fluid dynam- 
ics (CFD), the prediction of high-lift flow fields remains 
a challenge. The flow fields around multi-element air- 
foils are complex and are known to be dominated by 
different flow mechanisms at different operating condi- 
tions. 

A significant portion of the flow around multi- 
element airfoils is transitional over a wide range of op- 
erating condition. The onset location of the transition 
process varies with, for example, angles of attack and 
Reynolds numbers. In addition to boundary layer tran- 
sition, free shear flow transition is also likely to occur 
in, for example, the flow over the slat. 
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As the slat wake convects downstream, it may fur- 
ther interact with the transitional flow over the suc- 
ceeding elements, such as the flow over the main el- 
ement and the flap. As a result, an accurate predic- 
tion of the flow transition represents a major challenge 
for the computation of the flow over multi-element air- 
foils and the use of a transition model often becomes a 
discriminating factor in the outcomes of a prediction. 
Turbulence models are frequently used in the high-lift 
flow calculations. In this practice, turbulence models 
are “turned-on” at a designated location on the air- 
foils. The designated locations are usually correlated 
with the measured location of transition onset and the 
computational results can vary significantly with the 
assigned transition location. 

In this paper, the flow fields over a multi-element 
airfoil are calculated numerically using two two- 
equation turbulence models. 

These computational exercises are necessary in an 
effort to identify the limitation of turbulence models in 
the prediction of transitional flows. The computational 
results may also help to investigate areas of the flow 
where experimental data are difficult to obtain. 

The Navier-Stokes solver used in the study is de- 
scribed in the next section. The Chimera composite 
grids used are also presented. 

II Flow Solver 

The flow calculations presented here were performed 
with the INS2D code 4,5,6 . The code solves the two- 
dimensional, incompressible form of the Reynolds- 
averaged Navier Stokes equations using an artificial 
compressibility approach. The INS2D code is a fi- 
nite difference solver, using Roe’s third order upwind- 
biased, flux-difference splitting for the convective terms 
and a second -order central differencing for the viscous 
terms. The code is capable of solving both steady and 
unsteady problems with point-to-point matched grids 
or Chimera overset grids. The details of the INS2D 
code can be found in the references cited above. The 
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INS2D code has been applied extensively to the flow 
over various single- and multi-element airfoils 4,5,6 and 
has been shown to be a very reliable tool for such 
purposes. Generating the necessary grids for multi- 
element airfoil flow calculations can be very time- 
consuming regardless of the grid topology one chooses 
to use. In this study, overset grids were used and the 
grids were generated using the OVERMAGG 7 software. 
OVERMAGG is an automated script system specifi- 
cally designed to perform overset grid generation for 
multi-element airfoils. The use of OVERMAGG has 
resulted in a significant amount of saving in time dur- 
ing this study. 

III. Turbulence Models 


The SST model 8 and a k-e 9 ' 10 ' 11 model were used 
to study the flow field. The SST model is available as a 
model option in the INS2D code. The k — e model was 
recently implemented as a separate module to the code. 
There are two other one-equation model options in the 
INS2D code. The SST and the two one-equation mod- 
els have been applied to the flow over the 30P30N ge- 
ometry and the results have been reported previously 3 . 

The SST model was designed to take advantange of 
both the robustness of the k — lj type of models in the 
near wall region and the free-stream independence of 
the k — e type of models in the outer part of boundary 
layers and for free shear layers. The eddy-viscosity was 
formulated to account for the transport of the principal 
turbulent shear stress in boundary layers. The SST 
model has been shown to perform very well for many 
different flows. 

The k — e model used here represents a variant of 
the model developed by Shih and Lumley 9 , with the 
constant (^(=0.09) replaced by a variable C ^ formula- 
tion. The new formulation of C^, with an explicit de- 
pendence on the mean strain rate, was developed based 
on the realizability constraints of the Reynolds stresses. 
The model has been show n to work well for a wide range 
of flows 10,1 1 . 

The model equations for k and e are, 
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where 

c/ (,) = ^SijSij + tiiAj 

&ij — &ij 2 Cijk^k 

Qij is the mean rotation rate viewed in a rotating ref- 
erence frame with the angular velocity u>*. The param- 
eter A s is determined by 
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The damping function is defined by 

= [1. - exp(— (aiflfc + a 3 Rl + a 5 /?jj))]» (5) 


where 

a, = 1.7 x 10 -3 , a 3 =l(T 9 , a 5 = 5xl0 _1 ° 


The near-wall boundary conditions for the turbulent 
quantities are 

u 4 

k = 0.25u 2 and e = 0.251— (6) 

T v 

where u T denotes the surface friction velocity. 

IV. Geometry 

The McDonnell Douglas 30P-30N landing config- 
uration was used in this study. Figure 1 shows the 
geometry and the stations where profile data will be 
presented. Both the leading-edge slat and the trailing- 
edge flap have a deflection angle of 30°. The airfoil 
has been tested extensively at NASA Langley Low 
Turbulence Pressure Tunnel (LTPT) 1,12 . The geometry 
has been used as a test case in a CFD Challenge Work- 
shop held at NASA Langley in 1993 and many compu- 
tational results obtained by using different solvers and 
models have been reported 2,3,6 . 

The computational grid consists of seven zones and 
a total of about 103,000 grid points. Figure 2 shows 
the C grids around the slat, the main element, and the 
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flap, which consist of 175x45, 381 x 133, and 221x85 
points, respectively. Two H grids were used for the 
main-element cove and for the wake region of the flap. 
A background H grid extends to the top and the bot- 
tom walls of the wind tunnel. When the grid for one 
element overlaps another element(s), holes and outer 
boundaries are determined by using the PEGSUS code. 
The communicaiton between the overset grids are also 
handled by the PEGSUS code. 

Comprehensive grid-independence studies have 
been performed and reported previously using similar 
overset grids for the same geometry 3,6 . These and other 
calculations for the same geometry indicate that the 
grid point distributions in the wake regions can have 
certain effects on the velocity profiles in those regions 
and the wall pressure is not sensitive to the grid den- 
sity. These results also show that further refinement 
of the grid from the current level does not change the 
results appreciably. 

It should be noted that the turbulence models were 
applied without including any explicit account for the 
flow transition. 

Results will be presented for a Reynolds number of 
nine million and two angles of attack (AOA) of 8° and 
19°. 

V. Results and Discussions 

Figure 3 shows the comparison of the surface pres- 
sure coefficient distributions for the 19° case. Both 
models predict the surface pressure well with a slightly 
faster external flow predicted by the k — e model over 
the suction surface of the slat. 

For AOA=8°, the surface pressure coefficient distri- 
butions are shown in Figure 4. The overall agreement 
with the measurement is good for the main element and 
the flap. On the suction surface of the slat, the k — e 
model is predicting a faster external flow. 

Figure 5 shows the velocity profiles on the main el- 
ement and on the flap for A0A = 19°. At x/c — 0.1075, 
there is almost an uniform offset between the data and 
the predictions. This may have been caused by possible 
improper data calibration 2 . The predicted slat wake is 
located slightly further away from the main element 
surface at all three stations (r/c=0. 1075, 0.45, and 
0.85). As a result, the merging of the slat wake and the 
boundary layer is apparently yet to occur at x/c = 0.85. 
The boundary layer velocity profiles were predicted well 
by the k — e model, which resulted in a better predic- 
tion of the main wake immediatedly behind the main 
element trailing edge (Figure 5.d. x/c — 0.89817). 

On the flap (Figures 5.e and 5.f, x/c = 1.0321 and 


1.1125), there seems to be a better overall agreement 
between the k — e model and the measurement in terms 
of the development of the slat and the main-element 
wakes. None of the models gives a good prediction for 
the deficit of the main-element wake. 

Figure 6 shows the velocity profiles on the main el- 
ement for AOA=8°. The slat wake is predicted well 
at x/c — 0.1075. The predicted main element surface 
boundary layer is more developed than what the data 
has shown. At the later stations, x/c — 0.45 and 0.85, 
the predicted velocity profiles agree well the measure- 
ment. 

Figure 7 shows the surface skin friction coefficient 
distributions on the main element for AOA=19°. The 
computed results agree very well with the available 
measured data. 

It should be noted that correction terms for transi- 
tion modeling have been proposed for both models 2,3,13 . 
These corrections have not been included in the calcu- 
lations presented here. 

It is apparent that, based on the results shown here 
and those published in the previous studies, the two 
turbulence models used in this study can give a fairly 
accurate account of the outer inviscid flow field around 
the multi-element airfoil. The near wall viscous flow 
regions have not been predicted well. A measure com- 
monly adpoted to improve the prediction of the near 
wall viscous flow was to include some type of correction 
terms for the boundary layer transition, as the measure- 
ment did show a fairly large region of boundary layer 
transition on all three elements for a range of flow con- 
ditions. Various approaches have been proposed 2,3 and 
have shown a significant improvement. These results 
indicate that the slat flow prediciton plays a critical 
role in the overall flow computation. 

The flow field around the slat is complicated and 
difficult to measure due to the small sizes of the slat. 
The laminar/turbulent separation and reattchment in 
the lower pressure surface and the boundary layer tran- 
sition on the upper suction surface have been proposed 
as among the phenomena that need to be modeled cor- 
rectly to obtain reasonable predictions. Indeed, Rum- 
sey et al. 2 have shown that the velocity on the main el- 
ement can be predicted well if the calculated slat wake 
agrees with the measured wake to begin with. 

Since the Reynolds number is large, the initial slat 
wake development depends largely upon the boundary 
layer at the point of separation. 

In the following, we examine the boundary layer on 
the suction side of the slat by using the k — £ model 
results, hopeing to gain a better understanding of the 
possible physical mechanisms at play in the slat flow. 
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Figures 8 and 9 show the distributions of the veloc- 
ity mangitude and the turbulent kinetic energy profiles 
in the direction normal to the surface at three loca- 
tions on the suction side of the slat for A0A=8° and 
19°. Note that the velocity magnitude and the turbu- 
lent kinetic energy have been normalized by the local 
skin friction velocity. As the boundary layer on the 
slat surface develops, it interacts with the outer “invis- 
cid” flow. For AOA=19°, Figure 8.b shows that, as a 
result of this interaction, or mixing, of the low speed 
flow near the surface and the high speed flow in the 
outer stream, a region of flow with inflectional veloc- 
ity distribution appears. This mixing phenomana may 
be similar to that in the free shear layer. The mix- 
ing region is also detectible for AOA=8° in Figure 8. a, 
but becomes more apparent further downstream on the 
slat. 

Figure 9 shows that these inflectional velocity dis- 
tributions have produced local peaks of turbulent ki- 
netic energy in the “mixing” zone of the flow (x/c=- 
0.071 and 0.0053). For AOA=8°, the local peaks seem 
to be developing and, for AOA = 19°, their values have 
exceeded those in the buffer layer. Therefore, the vis- 
cous layer on the slat suction surface exhibits a double- 
deck structure, with a boundary layer near the wall and 
mixing-layer-like behavior in the outer region. 

As was mentioned earlier, there is no detailed flow 
profile measurement available for the slat region. The 
observation made here has been based entirely upon 
the computed results. However, to the extend that the 
observation is valid, it is evident that the transition in 
the shear layer and in the boundary layer are among 
the essential features of the flow around the slat. 
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Figure 9. Slat turbulent kinetic energy. (a) 8°; (b) 19°. 
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